A weighted adaptive least-squares finite element method for the Poisson–Boltzmann equation

The finite element methodology has become a standard framework for approximating the solution to the Poisson–Boltzmann equation in many biological applications. In this article, we examine the numerical efficacy of least-squares finite element methods for the linearized form of the equations. In particular, we highlight the utility of a first-order form, noting optimality, control of the flux variables, and flexibility in the formulation, including the choice of elements. We explore the impact of weighting and the choice of elements on conditioning and adaptive refinement. In a series of numerical experiments, we compare the finite element methods when applied to the problem of computing the solvation free energy for realistic molecules of varying size. The Poisson–Boltzmann equation (PBE) is a nonlinear, elliptic partial differential equation (PDE) used to model the elec-trostatic potential surrounding a fixed biomolecule immersed in an ionic solvent [1,2]. Calculation of the potential is an important ingredient in many molecular simulations, which necessitates fast, accurate numerical approximation of solutions to the PBE [3–5]. Yet, effective numerical methods for the PBE are challenging to construct due to the complex geometry, discontinuous coefficients, and singularities from point charges that are inherent in the model [6]. A wide variety of numerical schemes have been developed for approximating the solution to the PBE (e.g., see [6] for a survey). Here, we focus on the development and analysis of finite element strategies due to the well developed supporting variational theory. Traditionally the Rayleigh–Ritz or Galerkin finite element method has been the method-of-choice for either the linearized [7,8] or the non-linear Poisson–Boltzmann equation [9–11]. The Galerkin method is a viable approach, due its relative simplicity, strong theoretical base, and proven results in practice. Even so, there is an opportunity for more advanced approximations when considering other variational formulations of the problem. In this article, we consider a least-squares finite element method [12] for the linearized PBE due to its inexpensive and effective error estimation in adap-tive refinement, and because the first-order variables (or fluxes) are specifically expressed in the approximation; this is a useful attribute for many problems.